home *** CD-ROM | disk | FTP | other *** search
/ The Arsenal Files 4 / The Arsenal Files 4 (Arsenal Computer).ISO / ham / sattrk31.tgz / sattrack-3.1.tar / SatTrack / src / sattrack / satnute.c < prev    next >
C/C++ Source or Header  |  1995-03-16  |  20KB  |  478 lines
  1. /******************************************************************************/
  2. /*                                                                            */
  3. /*  Title       : satnute.c                                                   */
  4. /*  Author      : Manfred Bester                                              */
  5. /*  Date        : 17Aug94                                                     */
  6. /*  Last change : 15Mar95                                                     */
  7. /*                                                                            */
  8. /*  Synopsis    : This function block calculates the nutation effects.        */
  9. /*                                                                            */
  10. /*                                                                            */
  11. /*  SatTrack is Copyright (c) 1992, 1993, 1994, 1995 by Manfred Bester.       */
  12. /*  All Rights Reserved.                                                      */
  13. /*                                                                            */
  14. /*  Permission to use, copy, and distribute SatTrack and its documentation    */
  15. /*  in its entirety for educational, research and non-profit purposes,        */
  16. /*  without fee, and without a written agreement is hereby granted, provided  */
  17. /*  that the above copyright notice and the following three paragraphs appear */
  18. /*  in all copies. SatTrack may be modified for personal purposes, but        */
  19. /*  modified versions may NOT be distributed without prior consent of the     */
  20. /*  author.                                                                   */
  21. /*                                                                            */
  22. /*  Permission to incorporate this software into commercial products may be   */
  23. /*  obtained from the author, Dr. Manfred Bester, 1636 M. L. King Jr. Way,    */
  24. /*  Berkeley, CA 94709, USA. Note that distributing SatTrack 'bundled' in     */
  25. /*  with ANY product is considered to be a 'commercial purpose'.              */
  26. /*                                                                            */
  27. /*  IN NO EVENT SHALL THE AUTHOR BE LIABLE TO ANY PARTY FOR DIRECT, INDIRECT, */
  28. /*  SPECIAL, INCIDENTAL, OR CONSEQUENTIAL DAMAGES ARISING OUT OF THE USE OF   */
  29. /*  THIS SOFTWARE AND ITS DOCUMENTATION, EVEN IF THE AUTHOR HAS BEEN ADVISED  */
  30. /*  OF THE POSSIBILITY OF SUCH DAMAGE.                                        */
  31. /*                                                                            */
  32. /*  THE AUTHOR SPECIFICALLY DISCLAIMS ANY WARRANTIES, INCLUDING, BUT NOT      */
  33. /*  LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A   */
  34. /*  PARTICULAR PURPOSE. THE SOFTWARE PROVIDED HEREUNDER IS ON AN "AS IS"      */
  35. /*  BASIS, AND THE AUTHOR HAS NO OBLIGATIONS TO PROVIDE MAINTENANCE, SUPPORT, */
  36. /*  UPDATES, ENHANCEMENTS, OR MODIFICATIONS.                                  */
  37. /*                                                                            */
  38. /******************************************************************************/
  39.  
  40. #include <stdio.h>
  41. #include <math.h>
  42.  
  43. #ifndef STDLIB
  44. #include <stdlib.h>
  45. #endif
  46.  
  47. #include "satglobalsx.h"
  48. #include "sattrack.h"
  49.  
  50. /******************************************************************************/
  51. /*                                                                            */
  52. /* getNutation: calculates nutation in longitude and obliquity, and the       */
  53. /*              equation of the equinoxes                                     */
  54. /*                                                                            */
  55. /******************************************************************************/
  56.  
  57. void getNutation()
  58.  
  59. {
  60.     double tu, tusq, tucb, eps, trueEps, l, lp, f, d, o;
  61.  
  62. /******************************************************************************/
  63. /*                                                                            */
  64. /*  check if nutation quantities need to be updated                           */
  65. /*                                                                            */
  66. /******************************************************************************/
  67.  
  68.     if (fabs(julianDate - lastJulianDateNute) < NUTEUPDATEINT)
  69.         return;
  70.  
  71.     lastJulianDateNute = julianDate;
  72.  
  73.     tu   = (julianDate - JULDAT2000) / JULCENT;
  74.     tusq = tu*tu;
  75.     tucb = tu*tusq;
  76.  
  77. /******************************************************************************/
  78. /*                                                                            */
  79. /*  now calculate the mean obliquity of the ecliptic (eps), based on the      */
  80. /*  epoch J2000.0                                                             */
  81. /*                                                                            */
  82. /*  for reference see the "Astronomical Almanac", 1986, pages B24 - B31       */
  83. /*                                                                            */
  84. /*  the equation for calculation of eps has been taken from the               */
  85. /*  "Astronomical Almanac", 1984, page S21 (in the Supplement section)        */
  86. /*                                                                            */
  87. /******************************************************************************/
  88.  
  89.     eps = (84381.448 - 46.815*tu - 0.00059*tusq + 0.001813*tucb) * CASR;
  90.  
  91. /******************************************************************************/
  92. /*                                                                            */
  93. /*  now compute the nutation in longitude (totPsi) and in obliquity (totEps)  */
  94. /*  thus giving the true obliquity of the ecliptic (trueEps) and from that    */
  95. /*  the equation of the equinoxes, also based on the epoch J2000.0            */
  96. /*                                                                            */
  97. /*  calculation of the nutation parameters totPsi and totEps is performed     */
  98. /*  by using the "1980 IAU Theory of Nutation"                                */
  99. /*  (see "Astronomical Almanac", 1984, pages S21-S26)                         */
  100. /*                                                                            */
  101. /*  fundamental arguments of the Sun and the Moon are : l, lp, f, d, o        */
  102. /*  (time-dependent arguments)                                                */
  103. /*  (from the "Astronomical Almanac", 1984, page S26)                         */
  104. /*  the function reduce() reduces their values into the interval (0, 2pi)     */
  105. /*                                                                            */
  106. /*  the nutation in longitude (totPsi) and in obliquity (totEps) are updated  */
  107. /*  when the update time interval is elapsed                                  */
  108. /*                                                                            */
  109. /******************************************************************************/
  110.  
  111.     l   = ( 485866.733 + 1717915922.633*tu +31.310*tusq + 0.064*tucb) * CASR;
  112.     l   = reduce(l,ZERO,TWOPI);
  113.  
  114.     lp  = (1287099.804 +  129596581.224*tu - 0.577*tusq - 0.012*tucb) * CASR;
  115.     lp  = reduce(lp,ZERO,TWOPI);
  116.  
  117.     f   = ( 335778.877 + 1739527263.137*tu -13.257*tusq + 0.011*tucb) * CASR;
  118.     f   = reduce(f,ZERO,TWOPI);
  119.  
  120.     d   = (1072261.307 + 1602961601.328*tu - 6.891*tusq + 0.019*tucb) * CASR;
  121.     d   = reduce(d,ZERO,TWOPI);
  122.  
  123.     o   = ( 450160.280 -    6962890.539*tu + 7.455*tusq + 0.008*tucb) * CASR;
  124.     o   = reduce(o,ZERO,TWOPI);
  125.  
  126.     getNutationSeries(tu,l,lp,f,d,o);
  127.  
  128.     trueEps    = eps + totEps;                                       /* [rad] */
  129.     equEquinox = totPsi * cos(trueEps);                              /* [rad] */
  130.  
  131.     return; 
  132. }
  133.  
  134. /******************************************************************************/
  135. /*                                                                            */
  136. /* getNutationSeries: calculates sine and cosine series of nutation           */
  137. /*                                                                            */
  138. /******************************************************************************/
  139.  
  140. void getNutationSeries(tu,l,lp,f,d,o)
  141.  
  142. double tu, l, lp, f, d, o;
  143.  
  144. {
  145.     double b[106], sinb[106], c[106], cosc[106], bc[106];
  146.  
  147.     int n;
  148.  
  149. /******************************************************************************/
  150. /*                                                                            */
  151. /* coefficients and arguments of the sine and cosine terms for the Fourier    */
  152. /* series representation (106 terms) of the nutation in longitude (sine)      */
  153. /* and of the nutation in obliquity (cosine), respectively                    */
  154. /* (all known long and short-term variations are included)                    */
  155. /*                                                                            */
  156. /*                                                                            */
  157. /* totEps    = b(0)*cos(B(0)) + b(1)*cos(B(1)) + ... + b(105)*cos(B(105))     */
  158. /*                                                                            */
  159. /* with b(n) = aat(n,0) + aat(n,1)*tu  for n = 0,...,14   (time-dependent)    */
  160. /* and  B(n) = at(n,0)*l + at(n,1)*lp + at(n,2)*f + at(n,3)*d + at(n,4)*o     */
  161. /*                                                                            */
  162. /* with b(n) = aa(n,0)                 for n = 15,...,105 (time-independent)  */
  163. /* and  B(n) = a(n,0)*l + a(n,1)*lp + a(n,2)*f + a(n,3)*d + a(n,4)*o          */
  164. /*                                                                            */
  165. /*                                                                            */
  166. /* totPsi    = c(0)*sin(C(0)) + c(1)*sin(C(1)) + ... + c(105)*sin(C(105))     */
  167. /*                                                                            */
  168. /* with c(n) = aat(n,1) + aat(n,2)*tu  for n = 0,...,14   (time-dependent)    */
  169. /* and  c(n) = aa(n,1)                 for n = 15,...,105 (time-independent)  */
  170. /*                                                                            */
  171. /* and  C(n) = B(n)                                                           */
  172. /*                                                                            */
  173. /*                                                                            */
  174. /* first :  15 terms with time-dependent arguments                            */
  175. /*          (either nutation in longitude or in obliquity is time-dependent)  */
  176. /*                                                                            */
  177. /*          the data field at[][] contains the code for the linear combina-   */
  178. /*          tions of the arguments of Sun and Moon (l, lp, f, d, o)           */
  179. /*                                                                            */
  180. /******************************************************************************/
  181.  
  182.     static double at[15][5] = {
  183.          {  0,  0,  0,  0,  1},
  184.          {  0,  0,  0,  0,  2},
  185.          {  0,  0,  2, -2,  2},
  186.          {  0,  1,  0,  0,  0},
  187.          {  0,  1,  2, -2,  2},
  188.          {  0, -1,  2, -2,  2},
  189.          {  0,  0,  2, -2,  1},
  190.          {  0,  2,  0,  0,  0},
  191.          {  0,  2,  2, -2,  2},
  192.          {  0,  0,  2,  0,  2},
  193.          {  1,  0,  0,  0,  0},
  194.          {  0,  0,  2,  0,  1},
  195.          {  1,  0,  2,  0,  2},
  196.          {  1,  0,  0,  0,  1},
  197.          { -1,  0,  0,  0,  1}
  198.     };
  199.  
  200.     static double aat[15][4] = {
  201.          { -171996,  -174.2,   92025,     8.9},
  202.          {    2062,     0.2,    -895,     0.5},
  203.          {  -13187,    -1.6,    5736,    -3.1},
  204.          {    1426,    -3.4,      54,    -0.1},
  205.          {    -517,     1.2,     224,    -0.6},
  206.          {     217,    -0.5,     -95,     0.3},
  207.          {     129,     0.1,     -70,     0.0},
  208.          {      17,    -0.1,       0,     0.0},
  209.          {     -16,     0.1,       7,     0.0},
  210.          {   -2274,    -0.2,     977,    -0.5},
  211.          {     712,     0.1,      -7,     0.0},
  212.          {    -386,    -0.4,     200,     0.0},
  213.          {    -301,     0.0,     129,    -0.1},
  214.          {      63,     0.1,     -33,     0.0},
  215.          {     -58,    -0.1,      32,     0.0}
  216.     };
  217.  
  218. /******************************************************************************/
  219. /*                                                                            */
  220. /* second : 91 terms with non-time dependent arguments                        */
  221. /*                                                                            */
  222. /******************************************************************************/
  223.  
  224.     static double a[91][5] = {
  225.          { -2,  0,  2,  0,  1},
  226.          {  2,  0, -2,  0,  0},
  227.          { -2,  0,  2,  0,  2},
  228.          {  1, -1,  0, -1,  0},
  229.          {  0, -2,  2, -2,  1},
  230.          {  2,  0, -2,  0,  1},
  231.          {  2,  0,  0, -2,  0},
  232.          {  0,  0,  2, -2,  0},
  233.          {  0,  1,  0,  0,  1},
  234.          {  0, -1,  0,  0,  1},
  235.          { -2,  0,  0,  2,  1},
  236.          {  0, -1,  2, -2,  1},
  237.          {  2,  0,  0, -2,  1},
  238.          {  0,  1,  2, -2,  1},
  239.          {  1,  0,  0, -1,  0},
  240.          {  2,  1,  0, -2,  0},
  241.          {  0,  0, -2,  2,  1},
  242.          {  0,  1, -2,  2,  0},
  243.          {  0,  1,  0,  0,  2},
  244.          { -1,  0,  0,  1,  1},
  245.          {  0,  1,  2, -2,  0},
  246.          {  1,  0,  0, -2,  0},
  247.          { -1,  0,  2,  0,  2},
  248.          {  0,  0,  0,  2,  0},
  249.          { -1,  0,  2,  2,  2},
  250.          {  1,  0,  2,  0,  1},
  251.          {  0,  0,  2,  2,  2},
  252.          {  2,  0,  0,  0,  0},
  253.          {  1,  0,  2, -2,  2},
  254.          {  2,  0,  2,  0,  2},
  255.          {  0,  0,  2,  0,  0},
  256.          { -1,  0,  2,  0,  1},
  257.          { -1,  0,  0,  2,  1},
  258.          {  1,  0,  0, -2,  1},
  259.          { -1,  0,  2,  2,  1},
  260.          {  1,  1,  0, -2,  0},
  261.          {  0,  1,  2,  0,  2},
  262.          {  0, -1,  2,  0,  2},
  263.          {  1,  0,  2,  2,  2},
  264.          {  1,  0,  0,  2,  0},
  265.          {  2,  0,  2, -2,  2},
  266.          {  0,  0,  0,  2,  1},
  267.          {  0,  0,  2,  2,  1},
  268.          {  1,  0,  2, -2,  1},
  269.          {  0,  0,  0, -2,  1},
  270.          {  1, -1,  0,  0,  0},
  271.          {  2,  0,  2,  0,  1},
  272.          {  0,  1,  0, -2,  0},
  273.          {  1,  0, -2,  0,  0},
  274.          {  0,  0,  0,  1,  0},
  275.          {  1,  1,  0,  0,  0},
  276.          {  1,  0,  2,  0,  0},
  277.          {  1, -1,  2,  0,  2},
  278.          { -1, -1,  2,  2,  2},
  279.          { -2,  0,  0,  0,  1},
  280.          {  3,  0,  2,  0,  2},
  281.          {  0, -1,  2,  2,  2},
  282.          {  1,  1,  2,  0,  2},
  283.          { -1,  0,  2, -2,  1},
  284.          {  2,  0,  0,  0,  1},
  285.          {  1,  0,  0,  0,  2},
  286.          {  3,  0,  0,  0,  0},
  287.          {  0,  0,  2,  1,  2},
  288.          { -1,  0,  0,  0,  2},
  289.          {  1,  0,  0, -4,  0},
  290.          { -2,  0,  2,  2,  2},
  291.          { -1,  0,  2,  4,  2},
  292.          {  2,  0,  0, -4,  0},
  293.          {  1,  1,  2, -2,  2},
  294.          {  1,  0,  2,  2,  1},
  295.          { -2,  0,  2,  4,  2},
  296.          { -1,  0,  4,  0,  2},
  297.          {  1, -1,  0, -2,  0},
  298.          {  2,  0,  2, -2,  1},
  299.          {  2,  0,  2,  2,  2},
  300.          {  1,  0,  0,  2,  1},
  301.          {  0,  0,  4, -2,  2},
  302.          {  3,  0,  2, -2,  2},
  303.          {  1,  0,  2, -2,  0},
  304.          {  0,  1,  2,  0,  1},
  305.          { -1, -1,  0,  2,  1},
  306.          {  0,  0, -2,  0,  1},
  307.          {  0,  0,  2, -1,  2},
  308.          {  0,  1,  0,  2,  0},
  309.          {  1,  0, -2, -2,  0},
  310.          {  0, -1,  2,  0,  1},
  311.          {  1,  1,  0, -2,  1},
  312.          {  1,  0, -2,  2,  0},
  313.          {  2,  0,  0,  2,  0},
  314.          {  0,  0,  2,  4,  2},
  315.          {  0,  1,  0,  1,  0}
  316.     };
  317.  
  318.     static double aa[91][2] = {
  319.          {   46,   -24},
  320.          {   11,     0},
  321.          {   -3,     1},
  322.          {   -3,     0},
  323.          {   -2,     1},
  324.          {    1,     0},
  325.          {   48,     1},
  326.          {  -22,     0},
  327.          {  -15,     9},
  328.          {  -12,     6},
  329.          {   -6,     3},
  330.          {   -5,     3},
  331.          {    4,    -2},
  332.          {    4,    -2},
  333.          {   -4,     0},
  334.          {    1,     0},
  335.          {    1,     0},
  336.          {   -1,     0},
  337.          {    1,     0},
  338.          {    1,     0},
  339.          {   -1,     0},
  340.          { -158,    -1},
  341.          {  123,   -53},
  342.          {   63,    -2},
  343.          {  -59,    26},
  344.          {  -51,    27},
  345.          {  -38,    16},
  346.          {   29,    -1},
  347.          {   29,   -12},
  348.          {  -31,    13},
  349.          {   26,    -1},
  350.          {   21,   -10},
  351.          {   16,    -8},
  352.          {  -13,     7},
  353.          {  -10,     5},
  354.          {   -7,     0},
  355.          {    7,    -3},
  356.          {   -7,     3},
  357.          {   -8,     3},
  358.          {    6,     0},
  359.          {    6,    -3},
  360.          {   -6,     3},
  361.          {   -7,     3},
  362.          {    6,    -3},
  363.          {   -5,     3},
  364.          {    5,     0},
  365.          {   -5,     3},
  366.          {   -4,     0},
  367.          {    4,     0},
  368.          {   -4,     0},
  369.          {   -3,     0},
  370.          {    3,     0},
  371.          {   -3,     1},
  372.          {   -3,     1},
  373.          {   -2,     1},
  374.          {   -3,     1},
  375.          {   -3,     1},
  376.          {    2,    -1},
  377.          {   -2,     1},
  378.          {    2,    -1},
  379.          {   -2,     1},
  380.          {    2,     0},
  381.          {    2,    -1},
  382.          {    1,    -1},
  383.          {   -1,     0},
  384.          {    1,    -1},
  385.          {   -2,     1},
  386.          {   -1,     0},
  387.          {    1,    -1},
  388.          {   -1,     1},
  389.          {   -1,     1},
  390.          {    1,     0},
  391.          {    1,     0},
  392.          {    1,    -1},
  393.          {   -1,     0},
  394.          {   -1,     0},
  395.          {    1,     0},
  396.          {    1,     0},
  397.          {   -1,     0},
  398.          {    1,     0},
  399.          {    1,     0},
  400.          {   -1,     0},
  401.          {   -1,     0},
  402.          {   -1,     0},
  403.          {   -1,     0},
  404.          {   -1,     0},
  405.          {   -1,     0},
  406.          {   -1,     0},
  407.          {    1,     0},
  408.          {   -1,     0},
  409.          {    1,     0}
  410.     };
  411.  
  412. /******************************************************************************/
  413. /*                                                                            */
  414. /* calculate the linear combinations of the arguments of Sun and Moon         */
  415. /*                                                                            */
  416. /* B(n) = C(n) = bc(n) = ....                                                 */
  417. /*                                                                            */
  418. /******************************************************************************/
  419.  
  420.     for (n = 0; n <= 14; n++)
  421.     {
  422.         bc[n] = at[n][0]*l + at[n][1]*lp + at[n][2]*f + at[n][3]*d + at[n][4]*o;
  423.         bc[n] = reduce(bc[n],ZERO,TWOPI);
  424.     }
  425.  
  426.     for (n = 15; n <= 105; n++)
  427.     {
  428.         bc[n] = a[n-15][0]*l + a[n-15][1]*lp + a[n-15][2]*f
  429.                              + a[n-15][3]*d  + a[n-15][4]*o;
  430.         bc[n] = reduce(bc[n],ZERO,TWOPI);
  431.     }
  432.  
  433. /******************************************************************************/
  434. /*                                                                            */
  435. /* calculate the nutation in longitude (totPsi) (sine series) and             */
  436. /* the nutation in obliquity (totEps) (cosine series)                         */
  437. /*                                                                            */
  438. /******************************************************************************/
  439.  
  440.     for (n = 0; n <= 14; n++)
  441.     {
  442.         b[n]    = aat[n][0] + aat[n][1]*tu;
  443.         c[n]    = aat[n][2] + aat[n][3]*tu;
  444.     }
  445.  
  446.     for (n = 15; n <= 105; n++)
  447.     {
  448.         b[n]    = aa[n-15][0];
  449.         c[n]    = aa[n-15][1];
  450.     }
  451.  
  452.     for (n = 0; n <= 105; n++)
  453.     {
  454.         sinb[n] = sin (bc[n]);
  455.         cosc[n] = cos (bc[n]);
  456.     }
  457.  
  458.     totPsi = 0.0;
  459.     totEps = 0.0;
  460.  
  461.     for (n = 0; n <= 105; n++)
  462.     {
  463.         totPsi += b[n]*sinb[n];
  464.         totEps += c[n]*cosc[n];
  465.     }
  466.  
  467.     totPsi *= CASR / 10000.0;                                        /* [rad] */
  468.     totEps *= CASR / 10000.0;                                        /* [rad] */
  469.  
  470.     return;
  471. }
  472.  
  473. /******************************************************************************/
  474. /*                                                                            */
  475. /* End of function block satnute.c                                            */
  476. /*                                                                            */
  477. /******************************************************************************/
  478.